Analytic results and weighted Monte Carlo simulations for CDO pricing 
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We explore the possibilities of importance sampling in the Monte Carlo pricing of a structured 
credit derivative referred to as Collateralized Debt Obligation (CDO). Modeling a CDO contract 
is challenging, since it depends on a pool of (typically ~ 100) assets, Monte Carlo simulations 
are often the only feasible approach to pricing. Variance reduction techniques are therefore of 
great importance. This paper presents an exact analytic solution using Laplace-transform and MC 
importance sampling results for an easily tractable intensity-based model of the CDO, namely the 
compound Poissonian. Furthermore analytic formulae are derived for the reweighting efficiency. The 
computational gain is appealing, nevertheless, even in this basic scheme, a phase transition can be 
found, rendering some parameter regimes out of reach. A model-independent transform approach 
is also presented for CDO pricing. 



I. INTRODUCTION 

Econophysics literature, especially due to the 
availability of high-resolution stock exchange trad- 
ing data, has initially been concerned with interpret- 
ing equity stylized facts [1-3] and equity derivatives. 
The past decade however, has shown a tremendous 
rise in the trading volume of credit derivatives [4], 
i.e., products depending on an event like bankruptcy, 
default or changes in the credit rating of a company 
or government. The buyer of the protection against 
such an event transfers his credit risk to the seller, 
and pays a periodic fee in return, maximally un- 
til the maturity of the contract. Although in this 
setting, credit derivatives are instruments of risk re- 
duction, since it is not necessary to own, e.g., a bond 
of the companies of interest, they open ground for 
speculation, too. 

The simplest credit derivative is the Credit Default 
Swap (CDS), which is a swap transferring the risk 
of holding a fixed income product of a single com- 
pany, such as a bond. In case the company defaults 
on paying the bond coupons, the buyer of the CDS 
is entitled to the face value of the bond. The price 
(the periodic payment to the seller) of a CDS is re- 
ferred to as CDS spread. The higher the spread, the 
riskier the market deems investing in the company 
in question. 

The simplest credit derivative is the Credit Default 
Swap (CDS) , which is a swap transferring the risk of 
holding a fixed income product of a single company, 
such as a bond. In case the company defaults on 
paying the bond coupons, the buyer of the CDS is 
entitled to the face value of the bond. The price (the 



periodic payment to the seller) of a CDS is quoted in 
bps, i.e., 10""* of the nominal value of the contract, 
and is referred to as CDS spread. The higher the 
spread, the riskier the market deems investing in the 
company in question. Opposed to CDS-s, structured 
products depend on the status of many underlying 
assets (e.g., the bonds of many companies) which, 
due to the interwoven nature of business relation- 
ships and macroeconomic factors, have a complex 
correlation structure. The subprime mortgage crisis 
of 2007-2008, has shown that the rising volume of 
such contracts [4] can lead to unforeseen instabili- 
ties. 

The complexity of the stochastic models used for 
pricing structured products are often not analyti- 
cally solvable, Monte Carlo simulations are essen- 
tial tools for quantitative analysts. Much effort has 
been spent on improving different aspects of these 
MC simulations, like speed and accuracy. An addi- 
tional important task is calibrating the model pa- 
rameters to market-observable prices of benchmark 
instruments. One approach to the latter problem 
is reweighting MC paths gained using a prior prob- 
ability measure (weighted Monte Carlo, or WMC, 
Avellaneda et al. [5] and in the context of credit 
derivatives Cont et al. [6]). Our work also involves 
reweighting MC paths, but the goal is reducing the 
variance of the Monte Carlo estimates of the ex- 
pected cashflow, not calibration to market observ- 
ables. The reweighting scheme is based on the Ra- 
don -Nikodym derivative, and is also referred to as 
importance sampling (see Sec. 4.3 in [7]). 

In this paper, we are dealing with collateralized 
debt obligations, which are contingent on the default 
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status of the constituents of a reference portfolio, 
such as Markit iTraxx Europe or CDX NA IG[8]. 
The contract can basically be viewed as a combi- 
nation of many CDS-s, however, the net loss on 
the portfolio is cut up into smaller intervals termed 
tranches. The seller of a CDO tranche pays the 
excess loss on the portfolio above a threshold (at- 
tachment point of the tranche), up to a maximum 
value (detachment point), and receives in return a 
periodic payment from the buyer (proportional to 
the remaining width of the tranche), referred to as 
CDO tranche premium. Standard CDO tranches for 
the CDX NA IG series include the equity (0 - 3%), 
the mezzanine (3 — 7%, 7 — 15%) and the super se- 
nior (15—100%) tranches. The attachment points of 
the super senior tranches of other indices range from 
15% to 35%; in our work, we used the 30—100% shce 
as a representative super senior tranche. In the mar- 
ket, these products are quoted in either basis points 
(10~*, e.g., a spread of 100 basis points means the 
annual premium payment fraction, although pay- 
ments are typically made semi-annually), i.e., the 
value of the periodic payment, or, assuming a fixed 
premium, the value of the upfront payment (in % 
of the tranche notional). For the sake of simplicity, 
in this paper we assume zero upfront in each case 
considered. 

Note that the cost functions of buyer/seller are 
non-linear, thus, the dependence structure between 
portfolio elements plays a crucial role in pricing a 
CDO tranche. "Bottom-up" approaches, including 
the Gaussian copula model [9] which became infa- 
mous during the recent financial crisis [10], try to 
estimate this dependence structure, and price CDO's 
consistently with single name credit defaults swaps 
(CDS's, these depend on the default status of a sin- 
gle company). "Top-down" models, in contrast, deal 
directly with the aggregate loss on the portfolio, 
thereby decreasing the number of model parameters 
(in bottom-up approaches, this is done by introduc- 
ing homogeneity assumptions) and giving up infor- 
mation about component risks (but see the "random 
thinning" procedure in [11]). In this paper, we con- 
sider a simple top-down compound Poisson model in 
order to retain analytic tractability and demonstrate 
MC possibilities. 

The paper is organized as follows: Sections II A 
and II B summarize model details and the relevant 
quantities. Section III introduces a general method 
for CDO pricing for models including a constant in- 
terest rate and deduces analytic formulas for the 
cash flow of the CDO contract in the Poissonian 
case. Section IV turns to presenting the possibili- 
ties of the Monte Carlo simulation, and presents the 
path-reweighting technique, which, for this simple 
model, can again be analytically verified. 



II. BASIC CONCEPTS 

A. Collateralized Debt Obligation — The basic 
model 

Let us assume that the CDO is based on a pool of 
A'comp companies, which, for the sake of simplicity, 
corresponds to a total notional of A^ = 1 USD. 

This model is based on a single default process 
Df which is a compound Poisson process in the fol- 
lowing sense: the default events occur according to 
a simple Poisson process of intensity p, and during 
the i^^ event, a fraction Ji of the companies default. 
The jump sizes Ji are independent and identically 
distributed random variables of exponential distri- 
bution with parameter A, i.e., P {Ji < x) = l — e~^^. 
We also use the notation = A""'^ for the expected 
value of the jumps. We assume that { Ji}j are also 
independent from the jump times. 

We considered two natural ways to define the ac- 
tual loss process with values in [0, 1]. The first one 
is referred to as the linear specification given by 

Lf' :=min(A,l). 
The exponential specification 

Lr:=l-exp{-A} (1) 

is obtained by a smooth transformation from Dt. 
For small values of Dt, the two quantities and 
Dt are close, which is the typical case for the relevant 
parameter regimes. 

In what follows, we use the exponential specifica- 
tion (unless otherwise indicated) and denote the loss 
by Lt for simplicity, and we say that, at time t, 
an Lt proportion of the companies defaulted. (Note 
that in each specification, the portfolio loss is con- 
tinuous. This is a natural simplification for typical 
index portfolios where A^comp is 100.) 

The buyer of the a CDO tranche [a, d] makes peri- 
odic payments (called premium leg) proportional to 
the outstanding notional (the remaining width) on 
the tranche until either the maturity M (the expiry 
of the contract) is reached or the loss exceeds the 
detachment point. The seller of the protection pays 
the default leg after each default event, which is the 
increment of min{Lt, d} — min{Lf , a}. 

The default of a company does not mean that it 
becomes entirely worthless, a fraction f of its original 
value is recovered, i.e., the portfolio loss Lt increases 
by a 1 — f proportion of the jump which occurred at 
time t. In the simplest setting, the recovery f is a 
deterministic constant value f e [0, 1) (see [9, 12, 13] 
for example). In this paper, we simply assume that 
f = 0. 
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B. Relevant quantities 

Let the interest rate r be constant in time. For a 
tranche [a,c?], we denote by 



— min (Lf , d) 



min (Li, a) 



(2) 



the loss on this tranche at time t. The phrase tranche 
loss is often used in the literature for the total 
loss at the maturity. 

The default leg present value (def PV) of a tranche 
is approximately the expected present value of the 
tranche loss, more precisely, the increments of the 
loss are discounted. In mathematical terms, 



def PV : 




(3) 



which is meant as a Stieltjes integral. The depen- 
dence of the def PV on the tranche is suppressed in 
the notation. 

The premium leg present value (premPV) is the 
expected present value of the total amount of the 
periodic payment by the protection buyer. The an- 
nual payment is spread x ONt where the spread is 
given in basis points (bps) and fixed in the contract. 
ONt is the outstanding notional of the tranche [a, d\ 
at time t, i.e. 



ONt = d - 



which stand for the present value of the tranche loss 
and that of the total amount of premium leg paid by 
the protection buyer respectively. Note that these 
are random variables, and their expectations are 



E [X^'^'^^'ij = def PV, 
E (^X^P"^^")) = premPVlbp, 



compare with (3) and (4). 
We define 



a: = — ln(l — a), 
: = -ln(l --d). 



(7) 
(8) 



(9) 
(10) 



Hence, 



Due to (1), the event {Lt exceeds a} is equal to 
{Dt exceeds a}. The same hold with d and c) respec- 
tively. This notation serves to reduce the length of 
subsequent formulas. 



III. ANALYTIC APPROACH 

In the compound Poissonian case, we can derive 
explicit formulas for the relevant quantities defPV 
and premPVlbp. The expressions contain an in- 
finite series representation, which converges faster 
than exponential, therefore, our method provides a 
promising approach to CDO pricing. The basic idea 
of the computations is that we decompose the un- 
derlying expectation according to the first passage 
of certain levels of loss. 



premPV = spread x E 



^ONt dt 



(4) 



where the expectation on the right-hand side is de- 
noted by premPVlbp. 

The aim of CDO pricing is to give a good estimate 
of the fair value of the spread (and upfront) for 
given tranches. Therefore, the equation 

def PV = spread x premPVlbp (+upf ront) 

has to be satisfied, since the left-hand side is the 
expected income of the protection buyer, whereas 
the right-hand side is that of the protection seller. 
The problem is now reduced to finding the values 
defPV and premPVlbp. 

Throughout the paper, we will also use the nota- 
tion 



^(def) . ^ 



j^(prem) 



I'M 



e-'-'dC^ 



M 



e-''*ONt dt, 



(5) 
(6) 



A. The defPV and premPVlbp expressed with the 
first passage time 

The loss Lt is an exponential transformation of 
the compound Poisson process Dt, see (1), and the 
computations can be done in terms of Dt- Hence, 
we introduce 

Th := min{i > : A > /i} 

the first passage time of level h for Dt- It turns 
out that the quantities defPV and premPVlbp can 
be expressed by the integral of the function 

MKM) :=E(e-'-^'-l(Tft < M)) (11) 

which is a modified Laplace transform of the first 
passage time, and 1 (•) denotes the indicator func- 
tion. Therefore, knowing the function ipr is enough 
to determine def PV and premPVlbp and also the fair 
value of spread and upfront. 

It is not difficult to show that instead of the origi- 
nal definition of def PV in (3) where the increments of 
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loss are added, with a new approach, an integration 
along the vertical axis can be done, see also Figure 
1. We obtain 

defPV = E ^^rmin{t>^)■.Lt>x}^^ ^ ^^^^ 

= E (^^ e-''^'^l{Th < M) e-'' dh 



(12) 

after a change of variable under the integral sign 
{x — > h = — ln(l — a;)), which corresponds to consid- 
ering the compound Poisson process Dt itself instead 
of Lt. 
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Figure 1. A sample scenario of the loss process Lt- The 
jumps in the solid red line indicate the default events 
occurring according to a Poisson process of intensity p. 
Assuming continuous payment, the default leg present 
value def PV is the sum of the discounted increments of 
the loss process Lt (see Eq. (12)) and the the premium 
leg present value premPVlbp is related to the grey area 
(in case r = 0, they coincide, see Eq. (13)). (color online) 

Similarly, dehnition (4) of the premPVlbp is re- 
placed by 

premPVlbp — 

(I'd i'min{0<t<M:Lt>x} \ 
I i ' 

( ^0 ^min(Th,M) \ (13) 

= E / e-''dh 



forward manipulations, one obtains from (13) that 



premPVlbp = - (l - e"''^^) id - a) 

H 

r 

r 



ipo{h,M)e-''dh (14) 
ipr{h, M) e''' dh, 



with the unknown function ipr depending on the de- 
tails of the model. An important remark is that 
formulas (12) and (14) are valid for any distribution 
of the process Dt, not only for the compound Poisso- 
nian case, we assumed only a constant interest rate 
and continuous payment possibilities for both sides 
(it is simple to generalize the results for determinis- 
tic interest rate functions e^'* -> exp{— /J r (s) ds}, 
but we omit this possibility in the present paper). 
Using the present approach, for any distribution of 
Dt, it is enough to determine the function ipr for 
pricing a CDO. In the next subsection, we calculate 
iPr in the compound Poissonian case. 



B. Partial differential equation for the Laplace 
transform of the first passage time 

For the compound Poisson process Dt, a series 
representation of can be given as follows. In or- 
der to avoid later confusions, we fix the value of the 
interest rate r, and we suppress the subindex of (p. 
The expectation in (11) can be decomposed accord- 
ing to the time and the size of the first jump of the 
process Dt since these are independent exponential 
random variables with parameter p and A respec- 
tively. After change of variable, one can obtain the 
integral equation 



pe-(p+r)Afg-A/. 

I'M ( ,h \ 

XI e^^+^'MA/ e^''(p{x,y)dx + l] dy 











(15) 

using the memoryless property of the exponentials. 
Differentiating (15), we can deduce the partial dif- 
ferential equation 

dhMf + >^dM'P +{p + r)dh^ + Xrif = 0. (16) 



One boundary value is obviously 
(f{h,0) = 0. 



(17) 



In terms of Figme 1, formulas (4) and (13) give the 
indicated area in two different ways. After straight- 



One has to be more careful at the other one. By 
definition, the function is constant 1 along the line 
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h — 0, but the bivariate function ip{h, M) is not 
continuous here. Therefore, we redefine cp at the 
boundary, or, more precisely, we can say that the 
definition (11) is vahd only if /i > 0, and we extend 
the function continuously. Since as h X 0, the prob- 
ability that the first jump exceeds h tends to 1, the 
boundary value is 

</5(0, M) = lim ^{h, M) = / e-^'ype-pydy 

= (l - e-(^-)-) . (18) 

C. Solution of the PDE 

The equation (16) is a second order hyperbolic 
partial differential equation, which contains extra 
terms of lower order. One way of solving it is per- 
forming Laplace transformation in both variables. 
Let 

poo poo 

ipst-^ / e-''''e-''^ip{h,M)dhdM 
Jo Jo 

be the Laplace transform. We will also use the func- 
tions 



MM) 



e-"'ip{h,M) dh, 



ip{h,M) dM 



for computing the Laplace transforms of df^j^jip, dh^fi 
and OmV- They can be given by integration by 
parts. In the calculation, the Laplace transforms 
of the boundary values (17) and (18) also appear. 

The Laplace transform of the equation (16) is 
written as 



t{t + p + r) 



Stifst - t 



+ ip + r) [ sipst 



Xtips 



+ XrLfst = 0. 



t{t + p + r) 

The unknown function ipst can be expressed easily 
P 1 



'fist 



t st + Xt+ {p + r)3 + Xr' 



(19) 



The elimination of variable t can be done by using 
the identity 



,1-e- 



da; = 



1 

P{P + a) 



with a = {{p + r)s + Xr)/{s + A). We get 
P 1 



MM) = 



p + r s + 



\r 



p+r 

1 -exp [-{p + r)M + 



XpM\ 
s + x) 



(20) 

The second inversion is not completely obvious. The 
difficulty is that, in the second term in (20), the 
variable s appears in two different places: in the 
denominator of the prefactor 1 / (s + Xr / (p + r)) and 
in the exponential as well. 

We could use the general identity 



J f{y) dy^ da; 



J^e-'^P-^^y{x)dx 



(21) 



with /3 = —Xr/{p + r) to proceed. Then the problem 
reduces to finding the inverse Laplace transform of 
the function 



gs{M) = exp 



XpM 

'JiL 

p+r 



S + 



in the s variable where s occurs only once. It can 
be solved by considering the series expansion of the 
exponential and by performing the inversion for each 
term individually using the general formula 



{n-l)\ 



dx 



1 



{p + a)" ■ 



One may notice that for each term of the sum in 
the series expansion, we get functions of the form 
X H> cx^^^e~'^^. These are to be integrated by the 
left-hand side of (21) in place of /. Hence, we also 
use the following series representation of the lower 
incomplete gamma function: 



r-^e""* dt 



(n- 1)! 



1 



/ \k 



^ k\ 

k=0 



The resulting formula is 

Mh.M) = ^ . e-'^h-{p+r)M 
p + r 



oo 

E 

n=l 



{p + rY^Ar 'i-^ f Xph\^ 1 

k=0 



kV 
(22) 

It is not hard to see that the solution (22) indeed 
satisfies the equation (16) along with the boundary 
values (18) and (17). One more important special 
case is if r = 0. The (22) reduces to 



■pM 



k=0 



-\h 



{Xhl 

kl 



(23) 
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which can be verified intuitively as follows. The left- 
hand side of (23) is equal to P{Th < M) by defini- 
tion. The right-hand side is the sum of the weights of 
those trajectories of Dt which give rise to the event 
{T/i < M}. Assume that Dt has n jumps in the in- 
terval [0, Af]. The jump sizes can be generated by 
a Poisson point process with intensity A along the 
vertical axis. Th < M is satisfied if and only if this 
Poisson process has k point in [0, h] with < k < n. 
In the general r > case, the same factors as on 
the right-hand side of (23) appear in the formula, 
but there is no explicit probabilistic interpretation 
of (22). 

The benefit of the computations is that, with (22), 
the value of def PV and premPV arc known explicitly 
using (12) and (14). The expression (22) is com- 
putationally stable because of the factorial in the 
denominator. This new method gives an unbiased 
answer also for those tranches for which one can- 
not guarantee enough sample paths with the usual 
Monte Carlo simulation. In these parameter regimes 
(e.g., pricing a super senior tranche), the results can 
be compared to those of the weighted Monte Carlo 
simulation. 



IV. MONTE CARLO SIMULATION 

Fortunately, the simple compound Poisson model 
allows not only for a solution, but the reweight- 
ing can also be treated analytically. From the 
technical point of view, there are two approaches 
to extracting the expectations of interest from the 
computer-generated realizations of the model. The 
first, to which we refer to as "path-based", cal- 
culates each quantity of interest for all generated 
paths, and afterwards calculates the statistics of the 
gained datasets. The second, to which we refer to as 
"surface-based" , estimates the time-dependent prob- 
ability density function of the portfolio loss, i. e., a 
simple matrix, see Fig. 2 . This work deals with the 
former approach since it is more suitable for our pur- 
poses, which is getting a better estimate of the price 
of the more senior tranches, i. e., to reduce the num- 
ber of Monte Carlo paths necessary for obtaining a 
certain precision. 



A. Monte Carlo with paths 

The ith realization of the process is translated to 
financial values taking into account only the cash 




Time [years] 

Figure 2. Loss surface of the compound Poisson process 
with intensity p = 0.5 and exponential jump sizes of 
parameter A = 10 (expected value j =0.1) generated by 
10^ Monte Carlo paths. The Lt = stripe is not shown, 
because it contains much larger probabilities than the 
central region, and would therefore obstruct visibility, 
(color online) 

flows: 

xl'^'^J'£e~^'-{itf~eZ'.) (24) 
fe=i 

femax 

^(prem) ^ g-rt, Q^^^^ (25) 

fc=l 

with — ioi ^ii ■ • ■ I ifc„,ax ~ ^ representing an 
equally spaced time grid. The necessity of discretiza- 
tion is commonly a drawback of simulations, in our 
case, however, it is nearer to a real world scenario 
since CDO payments are typically transferred quar- 
terly. For efficiency reasons, the generated paths 
should be evaluated simultaneously for all tranches 
of interest. The individual financial values Xj*-^^*-* 

and X-^^^"^ per tranche are subject to a statistical 
analysis. 

B. Reweighted Monte Carlo 

The aim of reweighting is to reduce the computa- 
tional time needed for finding the expected present 
values of the premium and the default. Apply- 
ing reweighting, the program generates Monte Carlo 
paths with an altered parameter set of the compound 
Poisson model, that is, with an altered intensity p' 
and an altered expected jump size ^. Using pa- 
rameters which describe a relatively calm economic 
situation, senior tranches are typically not reached 
by MC paths, i. e., one obtains a poor estimate of 
their price. The idea is to simulate paths using a 
parameter set describing a more severe situation, 
and reweight them to preserve the original expected 
value while reducing the variance. 
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The reweighting relies on the Radon - Nikodym 
derivative, the derivative of one probabiHty measure 
with respect to another [7]. To get the appropri- 
ate weights, we have to calculate it for the real with 
respect to the altered one. 

We write the probability of realization of a path 
considering that the Nm jumps occur in small inter- 
vals {tj,tj -\- dt) j = 1, . . . ,Nm with the given in- 
tensity p and that there is no jump outside of these 
intervals. We handle the jump sizes in a similar way 
under the condition that the number of jumps is al- 
ready given by the Poisson process (we introduce the 
notations to = and Dq = 0): 



P(path) = Yl (e-''(*'-*'-i)pdi) • c'P'^^'-'^m) 

i=l 

Nm 

X n(e-^(^'.-^'.-iUd/i), (26) 
1=1 

with the first term corresponding to the pdf of the 
intervals between jumps, the second to the proba- 
bility that there is no event between the last jump 
and the maturity M and the third to the pdf of the 
jump sizes. It is straightforward to generalize the 
latter equation for arbitrary renewal processes, one 
only has to change the terms to the jump time and 
size distributions of interest. 

The Radon ~ Nikodym derivative of one probabil- 
ity measure with respect to another one is the ratio 
of the weights of the same path given in (26) under 
the two measures, and for a path with Nm jumps is 
given by 



R{Nm,Dm) 
_ dP 
~ dP 



-{Nm,Dm) 

Nm 



(27) 



,-ip-p')M-{\-\')DM 



where P and P' are respectively the measures 
parametrized by the real and the alternative param- 
eters. This quantity is both calculable analytically 
as a random variable and numerically for a specific 
Monte Carlo path. The Monte Carlo simulation cal- 
culates i?(path) step by step, during the generation 
of a path. Starting from the value 1, at each jump 
in the path, the program multiplies the stored value 
by the contribution of that jump (terms under the 
product signs in (26)); at the end of the path, it mul- 
tiplies the value by the contribution which describes 
that no more events happened until the maturity 
(middle factor on the right-hand side of (26)). 

In mathematical terms, the random variable sim- 
ulated under the altered measure P' is RX (with X 



standing for either premium or default leg), thus, its 
observable variance is 



Va.r'{RX) = ^'{R^X"^) - (E' RXf 
= E{RX^) - {EXf, 



(28) 



since E' {RX) = E {X), by definition of the Radon - 
Nikodym derivative. 

Being able to measure the variance, we use this 
feature to find the optimal parameter set for the 
speed of convergence and then we perform impor- 
tance sampling with those parameters. 



C. Reweighting: Analytic approach 

In this section, we analytically evaluate the vari- 
ance of the reweighted default considering the com- 
pound Poisson process Dt where the jump times 
follow a Poisson process with intensity p and the 
sizes of the jumps are independent exponentially dis- 
tributed random variables with parameter A, which 
are independent of the Poisson point process as well. 
Then we calculate the expected loss and its variance 
analytically with the assumption that there are no 
discount factors, i. e., r = 0. 

Recall (5), and note that in case r = 0, we have 

X^""^^ ^(Lt-a) • l(Lt e [a,d]) 
+ {d- a) ■ TL{Lt > d) 
= (1 -e-^' - a) • 1(A e [a,d]) 
+ (rf - a) • 1 (A > f ) • 



(29) 



Then 



defPV = E (^X^"^^^^^ = E' (^RX'-^"^' 



(30) 



where the second equality holds by the definition of 
the measure change. The important quantity here 
is the error of the Monte Carlo simulation carried 
out with importance sampling, thus, our aim is to 
calculate the variance given in (28). 

The joint density of Nm (number of jumps) and 
Dm is given by 

/* ("-, h) 

= P*(7Vm = n, Dm e ih,h + dh)) 

nl [n — 1)! 



n\{n — 1)! 



dh + o{dh) 



where * can stand for either altered or real. This 
bivariate joint probability density function is com- 
posed of the product of the probability density func- 
tions of a P0I(/9*M) describing n events until the 
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maturity and a r(A*,n) describing the conditional 
probability of arriving in {h, h + dh) having n jumps. 
Please note that this is a defective probability distri- 
bution, the missing mass is P*(Z?m = 0) = e^''**^. 
Using (29), one can calculate 



E = (1 - a) / Y.f{n,h)dh 

n=l 

e-''^/(n,/i)d/i 

n=l 



(31) 



+ {d~a) I ^/(n,/j)d/i, 



similarly, for the variance. 



3 oo 

J2e^^''R(.n,h)f{n,h) dh 

n=l 



-2(1 -a) / 

J a 

+ (d-a)2 



g oo 

^e-''i?(7i,;i)/(n,/.)d/. (32) 

n=l 

E 

n=l 



R{n,h)f{n,h) dh 



^R{n,h)f{n,h) dh 

n = l 



where the integrals can be expressed in terms of in- 
complete gamma functions, since the dependence of 
the integrands on /i is a product of a polynomial and 
an exponential function, that is, they are of the form 



-uh in 



dh. 



Deriving the result is not extremely difficult but 
rather technical, therefore we omit these details. 

We remark one more interesting feature of the ex- 
pectation in (32) which is the appearance of phase 
transition. In the first two terms on the right-hand 
side of (32), the integrals are not necessarily finite 
(for the other term, we do not have this issue). It 
can be verified by analyzing the exponential fac- 
tors in h of the integrand. In the first term on 
the right-hand side of (32), i?(n, h) contributes with 
g-(A-A )h ^ whereas f{n,h) gives an exponential fac- 
tor of e^^'^. The product of these two is clearly 
g-(2A-A jjgj^(,g^ fiYst integral in (32) is finite 
if and only if 



2A - A' > 



, 1 



but it does not mean a restriction in the practical 
point of view, because as mentioned earlier, one can 



only expect an improvement in the variance in case it 
is more likely to reach a senior tranche under P' than 
under P. (For the finiteness of the second integral, 
a weaker condition is sufficient.) 



D. Differences between the Monte Carlo and 
the analytic method 

We emphasize here that the difference between 
the analytic and the Monte Carlo methods is not 
merely technical (e.g. discretization), but conceptual 
as well. The analytic method gives all quantities of 
interest for a given input parameter set (if it corre- 
sponds to an analytically solvable model), namely, 
we get immediately expected value and variance for 
a given single tranche. 

In contrast, the Monte Carlo simulation has to 
generate a given number of simulations divided into 
packages (or not) to give the financial expected val- 
ues def PV and premPVlbp for a chosen tranche, fur- 
thermore, doing statistics on these experimental fi- 
nancial values provides information about their vari- 
ance. At this point the Monte Carlo method shows 
an advantage. Having all these realizations, we can 
simultaneously calculate the expected values and 
variances for any chosen set of tranches, without 
generating new paths. Hence, for efRcient measure- 
ment, one has to define a set of observed tranches 
and query the financial values simultaneously. 

As pointed out several times, besides the design 
differences, the Monte Carlo method provides a pow- 
erful tool in cases where no analytic solution can be 
found. Computers can still model and capture com- 
plex processes generating realizations, even if there 
are several specialties in the contract. 



E. Monte Carlo results 

In this section we discuss the simulation results 
and compare them to the analytic solution. For illus- 
tration purposes, we have chosen from the standard 
tranches (see Section I) the super senior tranche 
{a — 0.3, d = 1). For each Monte Carlo simulation 
we used 10^ paths and no recovery (f = 0). 

Approximating the model parameters in a calm 
economic situation by: 



0.05 



events 



year 



(33) 



for the intensity of the compound Poisson process 
and about 

1 ^ part of the original notional 
A ' event 
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for the expected number of defaults (i.e., 10 assets 
are expected to default per event for a 100-element 
portfolio). The standard maturity is 



M = 5 years. 



(35) 



The gain of the reweighting procedure is given by 
the well-known fact that the standard deviation of 
a sample average (of independent realizations of a 
random variable) is where a is the standard 
deviation of the random variable and n is the number 
of realizations used to estimate its mean-value. So, 
for the same precision, whatever we gain on a the 
square of it is gained on the number n of Monte 
Carlo paths to be generated 



n 
n' 



(36) 



The real-world parameters predict relatively rare 
events (one event forecast for four realizations) and 
small chance to touch the super senior tranche. 
Thus, one can expect the largest decrease in a for 
the super senior tranche when the alternative pa- 
rameters p' and ^ are chosen to be larger than the 
real ones which associates more trials to the super 
senior tranche without neglecting the others. 



Oh 
P3 



Q 



6.8 

6.6 

6.4 

6.2 
200 

100 





MC avg. of X 
E(X(''''*') istd.dev. 



0.0 0.2 0.4 0.6 0.8 1.0 



A L event J 



Figure 3. The measured average and standard deviation 
of X^^"^^ for the super senior tranche as a function of 
with f) = p, (A = 10, p = 0.05). The average and 
standard deviation of X^'''^'^ is given in basis points. The 
variance remarkably decreases for ■ 
(color onUne) 



2.8 X 



0.28. 



We executed a simulation for the given p, \ and 
M = 5 to study the effect of the jump size parameter 
on the standard deviation. The results for the super 
senior tranche are presented in Figures 3 and 4. The 
X(<i«*) and X^P""") scales are in basis points (BP), 
i.e., values are multiplied by 10"*, according to finan- 
cial convention. In the top part of each diagram we 
represent the empirical expected value, while in the 
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Figure 4. The measured average and standard deviation 
of X'P"") for the super senior tranche as a function of 
with p' — p, {\ = 10, p — 0.05). The average and 
standard deviation of X'^^'^'^"^ is given for a 1 basis point 
spread. The estimation of the premium leg present value 
was not improved for the super senior tranche. We have 
no analytic formula for premPVlbp. (color online) 
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Figure 5. The measured average and standard deviation 
super senior tranche as a function of p' , 
with A' = A, (A = 10, p — 0.05). The average and stan- 
dard deviation of the x'***^*' is given in basis points. The 
variance decreases most if p' G (4p, 5p) — (0.20,0.25). 
(color online) 



background the analytic one, together with the ana- 
lytically calculated Monte Carlo error as confidence 
interval. This error is given by the standard devia- 
tion of theX(^«*) for one path divided by ^/n, where 
n is the number of the Monte Carlo paths used in 
the measurement. The bottom part of the diagrams 
retraces the analytic and the empirical value of the 
same standard deviation to facilitate the compari- 
son. The expected gain in computational time is 
proportional to the decrease of the variance, that is, 
the square of the decrease in the standard deviation 
plotted here. 



10 



Pi 



a 



3.52 



3.50 



3.48 
10.0 

5.0 

0.0 



+ ' f » + +^-n-^+jf ++ ^ + ^ ++ 

MC avg. fo|: + 

E(X(P"°'')) rf 




0.0 



0.1 



0.2 



0.3 



0.4 



0.5 



/ events 
[ year J 



Figure 6. The measured average and standard deviation 
of for the super senior tranche as a function of 

p' , with A' = A, (A = 10, p = 0.05). The average and 
standard deviation of X'^"^™' is given for a 1 basis point 
spread. The estimation of the premium leg present value 
was not improved for the super senior tranche, (color 
online) 



We observe a very good match of the analytic and 
numerical (path) approach, and in addition, we dis- 
cover that we can realize gains in precision for the 
index tranche as well as for the tranches above the 
equity tranche. The figures confirm our expectation 
that the variance of X^^^^^^ can be decreased. The 
maximal gain for the super senior tranche is obtained 
for ^ — 0.28 where the variance was reduced to 
14% of its original value. In contrast, premPVlbp 
cannot be estimated more accurately than in the 
non-reweighted case. It should not discourage us, 
because the most uncertainty of the CDO contracts 
originates from the wrong estimation of X('^®*-''s 
properties. The properties of X^^^^'^^ can already be 
easily and accurately estimated without reweight- 
ing as shown in the figures. Even the reweighting- 
increased relative error ( ^^<^- dev. \ ^j? premPVlbp re- 

^ expected V. ^ ^ ^ 

mains an order below the relative error of def PV. 
The reweighting is therefore improving the most im- 
portant part of the pricing. 

A similar investigation was done on the jump time 
parameter (i.e., the intensity of the compound Pois- 
son process). Using the same real- world parameters 
the results for the super senior tranche are repre- 
sented in Figures 5 and 6. 

We observe an equally good match of the analytic 
and numerical (path) approach as previously. Here, 
we are glad to observe that all the tranches exhibit 
gain in precision with growing process intensity. The 
figures confirm our expectations that the variance of 
X^"^®*' can be decreased. The maximal gain for the 
index tranche is obtained for p' = 0.22 where the 
variance was reduced to 32% of its original value. 




0.0 1.0 2.0 3.0 4.0 5.0 6.0 

p' [events/year] 

Figure 7. The computational time as a function of the 
reweighting parameter p' , which denotes the Poisson 
process intensity. The points fit on a line of equation 
t — c + bp' where c = 1.1(67) is the average initialization 
cost of the 10*^ MC paths and b = 5.57(1) is the time 
consumed by processing the jumps, (color online) 



Although by increasing p we can gain on the num- 
ber of Monte Carlo paths, the gain on computation 
power is not so evident. On average, each path con- 
tains Mp' jumps which need to be generated by the 
random number generator and administrated by the 
financial layer. An explicit chronometry was done 
where the process' real computational time on the 
CPU was collected (not the clock ticks during the ex- 
ecution because the computer could do other things 
in the background). The time need for increased p' 
is plotted in Figure 7 and it is linear in p' . With re- 
gard to this fact the real gain in computational time 



Gtir 



tr 



comp 



(37) 



where ' stands for the reweighted simulation's re- 
sults, cr denotes the standard deviation of the mea- 
sured quantity [X'^^^^^ or X^p^®"^) and p denotes the 
Poisson process intensity. This effect is not present 
in the case of jj because higher jumps do not pro- 
voke more events - they only result in larger nu- 
meric values. We note here that the simulation with 
n — 10^ paths near the optimum region takes about 
3 seconds for each point. 

In possession of all this information, we expect 
that the optimum for speed will be somewhere in 
the region where both altered the intensity and the 
jump size parameter are higher than the real ones. 
To find the maximal gain, we have to analyze it in 
two dimensions, therefore we draw a map for the 
gain in MC paths and the gain in time as defined 
respectively in (36) and (37). 

If we look at Figures 8 and 9 we can conclude 
that the Monte Carlo reweighting method is success- 
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ful in sparing computational time. For the pricing 
of the super senior tranche the total saving in com- 
puter time reaches 91% of the original time without 
reweighting, which is more than appealing. (The 
total saving is defined as 100% • (1 — ^) where * 
stands for number of paths or time.) The saving 
in the number of paths was at least 66% for every 
tranche, worse for the equity tranche. The numerical 
results are summarized in Table I. 

With some further considerations, one can show 
that these figures should be smooth, which is satis- 
fied except for the region of rare but large jumps. 
This is explained by the relatively small number of 
paths (10^) being simulated, an order less than for 
Figures 3 to 6. This small number of 10^ paths is, 
however, justified, because the full map with this 
acceptable resolution is calculated in more than an 
hour on a state-of-the-art computer, the optimal 
sampling of this map is out of the scope of the cur- 
rent paper. We would like to assure the reader that 
the [0,0.5] margins on the maps were intentionally 
left out. We see already in the [0.5, 1] region the ten- 
dency of increasing variance and the analytical work 
suggests that the variance explodes in this model if 
the altered expected jump size goes below the half 
of the original. 

In addition to the above mentioned cases, we 
have analyzed some more extreme, crisis-like situ- 
ations where p ~ 1^^- It needed even more com- 
puter power because of the increase in intensity (thus 
jumps) in a M = 5 year CDO contract. As antici- 
pated, the gain in either number of paths or compu- 
tational time is less spectacular, since even the origi- 
nal estimation was not as poor as with the previous, 
"calm" parameter-set. This simulation showed that 
the reweighting has its limits, even in the number 
of paths less than 75% saving was achieved for the 
super senior tranche and at most 25% or nothing for 
the others, in contrast to the minimum of 66% of the 
previous case. 

V. CONCLUSIONS 

In this paper, we have shown the capabilities of 
importance sampling in estimating the fair price of 
CDO tranches. The simple model we covered, en- 
abled us to check our results both analytically and 
by computer simulation. We showed that this ap- 



proach is promising in pricing rare events, neverthe- 
less, it has to be treated with care, since even in this 
basic model, singular behavior emerged. 

Future directions include testing the method for 
more elaborate models (such as [14] or [15]) and au- 
tomatizing the optimization procedure, more specif- 
ically, gain a parameter set that simultaneously im- 




Figure 8. The magnitude of the gain — logj Gnum for the 
super senior tranche in the defPV measurement, given 
in number of Monte Carlo patlis as a function of the 
reweighting parameters: p denotes the altered Poisson 
process intensity, jj the altered expected jump size. The 
real process uses the real world parameters described in 
(33) to (35). The variance decreases for the super senior 
tranche by increasing the intensity and the jump size to 
2-6 times of their original value, (color online) 



proves all standard tranches. 
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The variance decreases for the super senior tranche by 
increasing the intensity and the jump size to 2-4 times 
of their original value, but this is less significant than in 
Figure 8. (color online) 
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